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We present a simple stochastic mechanism which generates pulse trains exhibiting a power law 
distribution of the pulse intervals and a 1//° power spectrum over several decades at low frequencies 
with a close to one. The essential ingredient of our model is a fluctuating threshold which performs 
a Brownian motion. Whenever an increasing potential V(t) hits the threshold, V(t) is reset to the 
origin and a pulse is emitted. We show that if V(t) increases linearly in time, the pulse intervals 
can be approximated by a random walk with multiplicative noise. Our model agrees with recent 
experiments in neurobiology and explains the high interpulse interval variability and the occurrence 
of 1//" noise observed in cortical neurons and earthquake data. 
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The omnipresence of 1//" noise in nature is one of the 
oldest puzzles in contemporary physics still lacking a gen- 
erally accepted explanation. The phenomenon is charac- 
terized by a certain behavior of the respective time signal: 
The power spectrum S{f) is proportional to 1//" at low 
frequencies / with a « 1. Examples include the light of 
quasars , electrical measurements |^ , music and speech 
1^, human cognition fl and coordination [0, the current 
through ion channels B , network traffic H] , burst errors 
in communication systems freeway traffic |9|, gran- 
ular flow [p^ , etc. The time signals of a large number 
of these systems P,[7|-[l^ resemble a pulse train consist- 
ing of individual, largely identical events which occur at 
discrete times. This is especially true for spike trains of 
single nerve cells for which 1//" noise has been observed 
in various brain structures The reported expo- 

nents, which depend both on the presence or absence of 
a sensory stimulus and on the state of the ani- 

mal (REM sleep vs awake state) [p"3|"p^, vary from 0.68 
to 1.38. The power-law behavior for spike train power 
spectra lies within the range 0.01 to 10 Hz, extending 
typically over 2 decades. In almost all cases the upper 
limit of the observed time over which fractal correlations 
exist is imposed by the duration of the recording. 

In this Letter, we propose a simple mechanism for gen- 
erating pulse trains with 1//" behavior in systems with 
a threshold-controlled dynamics like, e.g., neurons and 
earthquake faults. Our model is based on an integrate- 
and-fire (laF) mechanism and consists of a single unit 
characterized by two variables (see Fig. ||): The voltage 
V{t) and the threshold C{t). Initially, the voltage is be- 
low the threshold. Then, the voltage increases monoton- 
ically in time — in the simplest case just linearly — and 
the threshold evolves according to a Brownian motion 
with diffusion constant D within reflecting boundaries 
Vq < Ci < C{t) < Cu- As soon as V[t) has reached the 
threshold, the voltage is reset to Vq and a pulse of unit 
height is emitted. In this way, a pulse train is generated. 
Note that the threshold is not reset to its initial value. 



Such a model is often used to describe single neurons: 
V{t) is the membrane potential and the emitted pulse 
is the generated action potential. However, it is usually 
assumed that the threshold is constant in time. Recent 
investigation have shown that this is not true for cortical 
neurons in vivo |]l6| and in vitro [IT^ . Additionally, there 
is evidence that the spike trains of auditory neurons |l^ ] 
and of neurons in the mesencephalic reticular formation 
|l9| are not renewal, i.e., successive time intervals be- 
tween spikes are correlated. These facts are incorporated 
in our model in the simplest possible way. Moreover, the 
effects of dead time or absolute refractoriness, which lim- 
its the rate at which a neuron can fire, are automatically 
included in our model via the lower bound C/. The up- 
per bound Cu prevents an infinite time difference between 
two pulses. 
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FIG. 1. Dynamics of our model for a linear increasing volt- 
age V{t) with Vb = 0. The dashed lines represent the lower 
and the upper boundary for the fluctuating threshold C{t). 

As one can see from Fig. the power spectrum of 
the pulse train generated by our model with linear in- 
creasing voltage shows a 1//" decay over several decades 
pO|. The frequency below which white noise behavior 
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is observed is determined by C„ — Ci and goes to zero 
for Cu — Ci — > oo. The exponent a is not universal and 
increases if the ratio (C/ — V^jy/D increases. In the 
limit {Ci - Vo) < a/D, we find a = 0.5. This result is 
explained by the fact that this limit corresponds to the 
case in which the waiting times between pulses are purely 
determined by the diffusive dynamics of the threshold. 
Hence, each waiting time is merely the first return time 
of a Brownian motion which obeys a power-law distribu- 
tion with exponent —1.5. This implies a = 0.5 [ pT[ . 

To compare our results with real neurons, consider the 
parameters of curve (b): The dead time of a neuron is 
typically of order of milliseconds and the maximal time 
difference between two spikes of the order of seconds. 
This is exactly the ratio between C; and C„. Moreover, 
we can now identify one unit of time in our model with 5 
milliseconds real time. Hence, the 1//" behavior in (b) 
is found for / < 11 Hz. This and the fact that a ~ 1.02 
reproduces the experimental results very well. 

Extensive numerical investigations have shown that 
our model is very stable with respect to variations of 
the dynamics. Different forms of voltage increase, e.g., 
a linear, a squared, or a square-root increase, give sim- 
ilar results. The assumption of a monotonical increase 
of the voltage can also be dropped. A small amount of 
noise can be added to the voltage signal without altering 
our findings. Substituting the reflecting boundaries by 
a confining potential does not change our results, cither. 
This points towards a generic behavior. 

To obtain the power spectrum of our model analyti- 
cally, consider the signal X(t) ~ 8[t — tk) where 
denotes the time of occurrence of the fcth pulse. It follows 
as shown in lE^ 
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With / = limT- 
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where (• • •) denotes the average over the ensemble and 
over k. Hence, we need to know the probability distri- 
butions 'i/q{T)dT of the time differences between pulses 
T = tk — tk+q for all integers q. V&g is merely the gth 
passage times density function gq averaged over the sta- 
tionary probability distribution of "initial" states V{0) 
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Since gq can be computed from the first passage times 
density function (FPTDF), we will first focus on the lat- 
ter. 
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FIG. 2. Power spectrum of the pulse train generated by 
our model. Parameters are (a): Vb = 0, Ci — 0.2, Cu = 40, 
\/0 = 0.2. (b) as (a) with C„ = 200. (c): Vb = 0, Ci = 0.02, 
Cu = 4000, VD = 2.0. This curve is shifted down by 1 
decade. All curves show a clear 1//" decay, a varies from 0.6 
to 1.1. 

The FPDTF of our model with linear voltage increase 
can be obtained by mapping the model to an laF model 
with constant threshold C = {Cu + Ci)/2: The voltage 
V{t) is defined as the sum of V{t) and C - C(t). This 
means that V{t) fluctuates around V{t). It also implies 
that V{t) is reset to Vb + C— C(tfc) after the fcth threshold 
crossing. Hence, the correlations are now encoded in the 
fluctuating reset. In conclusion, V{t) behaves almost as 
a Brownian motion with drift to an absorbing barrier 
with a reset following each barrier crossing. The only 
difference is that the stochastic process is restricted to 
the interval [V{t) + C - Cu,V{t) + C - C'l] at time t 
which makes an exact mathematical treatment difficult. 
Neglecting the restriction for a moment and setting C = 
0, we obtain from the associated Fokker-Planck equation 
the FPTDF 
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and gi = for t < 0. Here, 1^(0) is the initial dis- 
tance. Eq. (j^) is, of course, just an approximation to the 
FPTDF of our model. This can already be seen from the 
fact that the FPTDF of our model is identically zero for 
T < (C; — Vq) and r > (C„ — Vq). However, as will be 
shown below this approximation proves to be very useful 
and Eq. (|) can even be simplified to a Gaussian with 
the same mean and variance: 
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From the FPTDF one could now compute higher passage 
times density functions gn in principle by convolution 

5„(t|1/(0)) = / gi{t\V{0))gn^i{T-t\t)dt. (7) 
Jo 
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However, this convolution resists an analytical treatment 
due to the Markovian character of our process and the 
form of Eqs. (^,0). The former manifests itself in the 
fact that the first passage time directly gives the initial 
distance between the threshold and the voltage for the 
next passage problem. The reason for this is that the 
norm of the kth reset of the voltage V equals the time 
difference between the fcth and the {k — l)th pulse due to 
the linear increase of the voltage V with slope 1. Hence 
the reset depends on the last passage time only and the 
whole stochastic point process is totally described by the 
FPTDF plus its Markovian property. 
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FIG. 3. Power spectrum of the time signal generated by the 
linear version of our model and by the two approximations for 
the same parameters as in Fig. ^, (b). 

These properties enable us to interpret the stochastic 
process generated by the Gaussian approximation (^) as 
a random walk of the inter-spike- intervals (ISI) : 



Tk+l 



Tk 



(8) 



Here, rjk denotes the white noise source. Note the spe- 
cial kind of multiplicative noise which distinguishes our 
model from the one in The variance of the step 

length is proportional to the current "position" . Hence, 
the origin is a fixed point of the random walk. To omit 
this difficulty and in spirit of our original model, we con- 
sider the random walk Tk to be confined by two reflecting 
boundaries, i.e., < (C„ - Vq) < Tk < (Ci - Vq). The 
power spectrum of the pulse-train generated by such a 
random walk with tk+i = tk + Tk shows a clear 1//" 
behavior with the same exponent as for our model (see 
Fig. H). This is also true if the new inter-spike-interval 
is chosen from an inverse Gaussian distribution given in 
Eq. with an initial condition depending on the last 
interval. Hence, these approximations seem to be jus- 
tified. This is further confirmed by the stationary ISI 
distribution. For the random walk, the ISI distribution 
function P{t) is proportional to t~^ |2^]. Simulations 
show that this is in excellent agreement with our model 
and the inverse Gaussian approximation. 



We have to point out that the behavior of P{t) de- 
pends on the specific kind of voltage increase. For a volt- 
age increase with {t — tiast)^ with < /3 < 2, we find 
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for our model. Substituting the reflect- 



ing boundaries by a potential of the form a/x -\- bx^ for 
example, we still find power law tails in the ISI distribu- 
tion. This is exactly what measurements show for corti- 
cal neurons pq] , especially in the mesencephalic reticular 
formation [|19[. Hence, our model seems to be especially 
well-suited to explain the occurrence of 1//" noise in 
that formation. In general, laF models can not explain 
the high interspike interval variability exhibited by corti- 
cal neurons |^ . A fluctuating threshold as described by 
our model, however, can solve this long-standing issue. 

For renewal processes the ISI distribution function and 
the FPTDF are identical and completely describe the 
process. This is the case for a standard laF neuron with 
reset of the voltage to the origin after each barrier cross- 
ing which is described by the FPTDF in Eq. (||). Such a 
model can neither explain a 1//" signal nor a power-law 
decay of the ISI distribution. In contrast to our model. 
Usher and co-workers showed that fractal behavior might 
be a consequence of the global activity dynamics of a 
network of laF neurons . Due to experimental limita- 
tions, however, the link between 1//" single-unit power 
spectra and macroscopic activity dynamics remains, as of 
now, a conjecture. Another attempt to explain the phe- 
nomenon of 1//" noise is based on fractal and fractal-rate 
stochastic point processes Certain types of these 

processes properly characterize the statistical properties 
found in different experiments. However, in many cases it 
is not clear a priori why the real system should generate 
such a process. This is especially true for the clustering 
Poisson process which was applied to explain 1 //" noise 
in the mesencephalic reticular formation ||l^ . In con- 
trast to that, our model provides a simple explanation of 
1//" noise in integrate-and-fire systems in general and, 
hence, applies to neurons in particular. The crucial as- 
sumption is a fluctuating threshold. Such a behavior is 
related to models presented in |29| and was already con- 
sidered in [^5|] to explain the power law decay of neuronal 
ISI distribution functions. However, the stochastic point 
process was still assumed to be renewal and could not 
explain a 1//" behavior of the power spectral density 
function. Consequently, the second crucial assumption 
of our model is the Markovian character of the process 
in agreement with Refs. [^^. To directly verify our 
model for neurons in the mesencephalic reticular forma- 
tion, one should measure the evolution of the reset from 
the threshold potential. 

All three main characteristics of our model, i.e. , a 1//" 
spectrum, a power-law decay of the ISI distribution, and 
a Markovian dependence of the ISI, can not only be ob- 
served in nerve cells but also in other systems. For ex- 
ample, the pulse signal, defined by associating the oc- 
currence of pulses of unit height with earthquakes in the 
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Mojave region |30|, shows 1//" fluctuations with a = 1.3 
and has an ISI distribution with exponent —1.1 (Fig. ||). 
Moreover, there is evidence that the ISI do not obey a 
renewal process. Rather a Markovian dynamics can be 
found Q. These findings imply that the effective evo- 
lution of the ISI is similar to the one generated by our 
model and caricatured by Eq. (g). 




f/Hz 

FIG. 4. Power spectrum of the earthquake signal. Inset: 
ISI distribution. 

To summarize, we have presented a simple stochastic 
model which is able to transform integrated white noise 
with a 1//^ power spectrum into noise with 1// tail. 
The basic mechanism is similar to what is expected to 
describe the dynamics of single neurons and, thus, our 
model can explain the behavior of neurons in the central- 
nervous-system 1 14 19]. The analysis of earthquake data 



even suggests that the effective dynamics of the ISI is the 
same for many systems exhibiting 1//" noise. 
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